function ROIout = findROI(path, imname, roisz, xOffset)
ROIout = zeros(roisz,roisz);
A = imread([path imname]);

G = fspecial('gaussian',[5 5], 2);      % Applying the Smoothing Filter
C = imfilter(A, G, 'same');
[w, h] = size(C);


BW = im2bw(C, graythresh(C));           % Convert output into Binary
[B,~,~] = bwboundaries(BW);             % Detect the Boundary

centroid = regionprops(BW,'Centroid');  % Locate the centroid
centroid = centroid.Centroid;
% Cendroid(1) - X coordinate of Centroid
% Cendroid(2) - Y coordinate of Centroid

outline = flipud(B{1});     % This is the edge we are interested in

regionLine = zeros(length(outline),3);  % Put all coordinates and distances in the same array

for i = 1:length(outline)
    regionLine(i,:) = [outline(i,2) outline(i,1) sqrt((outline(i,2)-centroid(1))^2+(outline(i,1)-centroid(2))^2)];
end

MINPEAKDISTANCE=10;

if length(regionLine)<MINPEAKDISTANCE*3
    return;
end
% Use -regionLine (flip the line up side down) to find minima...
[~, i] = findpeaks(smooth(-regionLine(:,3),50),'MINPEAKDISTANCE',MINPEAKDISTANCE,'MINPEAKHEIGHT',mean(-regionLine(:,3)));
% Higher order Moving average filter used... User may wish to change the
% parameters of the above function

if length(i)<4
    return;
end

if numel(regexp(imname,'_l_')) %Left Image minima of interest 1:2 and end-1:end
    coordinates = [regionLine(i(1:2),1:2); regionLine(i(end-1:end),1:2)];
else %Right Image minima of interest 1:3 and end
    coordinates = [regionLine(i(1:3),1:2); regionLine(i(end),1:2)];
end

sortCoord = sortrows(coordinates,2);    % Sort the coordinates from top of image to bottom

sc = lineScore(sortCoord);              % Distance based score to reject one point
if sc(1)<sc(2)
    sortCoord(4,:) = [];                % Reject point 4 (1, 2, 3 are close)
else
    sortCoord(1,:) = [];                % Reject point 1 (2, 3, 4 are close)
end

% Keeping track of the centroid after rotation
zeroIm1 = zeros(w,h);
zeroIm2 = zeros(w,h);
% zeroIm3 = zeros(w,h);

if length(sortCoord)<3
    return;
end

radius=2;

if round(sortCoord(1,2))<radius+1
    sortCoord(1,2)=radius+1;
end
if round(sortCoord(1,1))<radius+1
    sortCoord(1,1)=radius+1;
end
if round(sortCoord(3,2))<radius+1
    sortCoord(3,2)=radius+1;
end
if round(sortCoord(3,1))<radius+1
    sortCoord(3,1)=radius+1;
end

if round(sortCoord(1,2))>w-radius
    sortCoord(1,2)=w-radius;
end
if round(sortCoord(1,1))>h-radius
    sortCoord(1,1)=h-radius;
end
if round(sortCoord(3,2))>w-radius
    sortCoord(3,2)=w-radius;
end
if round(sortCoord(3,1))>h-radius
    sortCoord(3,1)=h-radius;
end

zeroIm1(round(sortCoord(1,2))-radius:round(sortCoord(1,2))+radius,round(sortCoord(1,1))-radius:round(sortCoord(1,1))+radius) = 255;
zeroIm2(round(sortCoord(3,2))-radius:round(sortCoord(3,2))+radius,round(sortCoord(3,1))-radius:round(sortCoord(3,1))+radius) = 255;
% zeroIm3(centroid(1)-2:centroid(1)+2,centroid(2)-2:centroid(2)+2) = 255;

% Angle of rotation to keep point 1 and 3 in the horizontal...
grad = (atan((sortCoord(3,2)-sortCoord(1,2))/(sortCoord(3,1)-sortCoord(1,1))))*180/pi;
% grad = 0;
grad_offset = 90; % 90
assignin('base','angle',grad);
if numel(regexp(imname,'_l_')) %Left Image
    rotPalm = imrotate(A,grad+grad_offset);          % Rotation Angle
    rotzeroIm1 = imrotate(zeroIm1,grad+grad_offset);
    rotzeroIm2 = imrotate(zeroIm2,grad+grad_offset);
    %     rotzeroIm3 = imrotate(zeroIm3,grad+90);
else
    rotPalm = imrotate(A,grad-grad_offset);          % Rotation Angle
    rotzeroIm1 = imrotate(zeroIm1,grad-grad_offset);
    rotzeroIm2 = imrotate(zeroIm2,grad-grad_offset);
    %     rotzeroIm3 = imrotate(zeroIm3,grad-90);
end

[top2, top1] = find(rotzeroIm1==255);          % Find Top coordinate in Rotated Image
[bottom2, bottom1] = find(rotzeroIm2==255);       % Find Bottom coordinate in Rotated Image

new_c = [round((top1(1)+bottom1(1))/2)+xOffset round(((min(top2)+max(top2))/2+bottom2(1))/2)];
% new_c = [xOffset(1) yOffset(1)];

% assignin('base','new_c',new_c)
% assignin('base','rot_p',rotPalm)

try
    ROIout = rotPalm(new_c(2)-round(roisz/2)+1:new_c(2)+round(roisz/2),new_c(1)-round(roisz/2)+1:new_c(1)+round(roisz/2));
    
catch err
    if (strcmp(err.identifier,'MATLAB:badsubscript'))
        return
    else
        rethrow(err)
    end
end
%% Showing an output - Comment/Uncomment as needed
% subplot(2,2,1:2)
% imshow(A)
% title(imname)
% 
% 
% subplot 223
% imshow(rotPalm)
% hold on
% plot(new_c(1),new_c(2),'b*')
% rectangle('Position',[new_c(1)-round(roisz/2),new_c(2)-round(roisz/2),roisz,roisz])
% plot(top1(1),top2(1),'r*')
% plot(bottom1(1),bottom2(1),'g*')
% title('Normalized and Rotated')
% 
% subplot 224
% imshow(ROIout)
% title('Region of Interest')
% drawnow
